Stability of the hard-sphere icosahedral quasilattice 



H. M. Cataldo 
Departamento de Fisica 
Facultad de Ciencias Exactas y Naturales 
Universidad de Buenos Aires 
14-28 Buenos Aires, Argentina 

C. F. Tejero 
Facultad de Ciencias Fisicas 
Universidad Complutense de Madrid 
E- 28040 Madrid, Spain 

Abstract 

The stability of the hard-sphere icosahedral quasilattice is analyzed using the 
differential formulation of the generalized effective liquid approximation. We 
find that the icosahedral quasilattice is metastable with respect to the hard- 
sphere crystal structures. Our results agree with recent findings by McCarley 
and Ashcroft [Phys. Rev. B 49, 15600 (1994)] carried out using the modified 
weighted density approximation. 
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I. INTRODUCTION 



The observation in 1984 by Shechtman et al! of a sharp diffraction pattern in a AIMn 
alloy with the symmetry of the icosahedron, opened a new field in condensed matter physics. 
Ever since experimental evidence of other materials having sharp diffraction patterns with 
symmetries forbidden by classical crystallography has continued to grow. The name qua- 
sicrystals has been coined to represent systems with perfect order but without periodicity, 
i.e. quasiperiodic systems. 

Since the nonperiodic three-dimensional (3-D) Penrose tiling! has a diffraction pattern 
closely similar to that of icosahedral alloys, it has been extensively studied to account for 
icosahedral point symmetry and also because of its relative simplicity. The 3D Penrose 
tiling is usually constructed by projection from a 6D simple cubic lattice. The projection is 
performed by first defining an acceptance domain in the 3D complementary space in order 
to select what points of the 6D simple cubic lattice are effectively projected to form the 
3D quasilattice. In the present investigation the 3D Penrose tiling has been generated by a 
special choice of the shape and size of the acceptance domain as described by Elseri. The 
last step to model a quasicrystal concerns the decoration i.e., the location of the lattice 
points forming the quasicrystal and a choice for the pair interaction potential. 

Of the many questions about quasicrystals, one concerns the stability of these phases. 
The first theoretical approaches to such a question!'! were based upon the Landau theory 
of crystallization where the free energy is expanded in powers of an order parameter related 
to density waves with icosahedral symmetry. As a main result it was shown that multi- 
component systems are required to achieve stability, a fact that agrees with the experimental 
findings!. 

On the other hand and from a more microscopic viewpoint, the analysis of the stability 
of quasicrystals may be rather difficult on general grounds, but an important simplifica- 
tion occurs if a hard-sphere pair potential is assumed. Indeed, a simple calculation of the 
maximum packing fraction, i.e. the fraction of the total volume occupied by the spheres, 
provides an important criterion of the stability chances. For instance, a 3D Penrose tiling 
with all vertices occupied by identical hard spheres leads to a fluid-like packing fraction and 
thus it must certainly be discarded as a model of a quasicrystal. Moreover, the interest of 
considering hard-sphere quasicrystals goes far beyond simplicity. Thus, numerical studied! 
have shown that the crucial criterion for the quasicrystal stability with more realistic, e.g., 
Lennard- Jones, interactions is the packing fraction of the quasicrystalline hard-sphere dec- 
oration. 

There are basically two ways of improving the poor packing fraction of the above fully 
occupied Penrose tiling. These options are to change the decoration! or the acceptance 
domain!. Both procedures give approximately the same optimal packing fraction ~ 0.63, a 
value that now indeed justifies a further stability study. Recently, McCarley and Ashcrofti 
have studied the hard-sphere quasicrystal using a modification of the acceptance domain!, 
to obtain from the modified weighted density approximation a metastable quasicrystal with 
respect to the crystalline and fluid phases. Their method is entirely formulated in the 6D 
reciprocal space which avoids direct summations over the quasilattice. However, a drawback 
of this method is that a truncation of the sum in the 6D reciprocal lattice is needed, the 
induced estimate error in the free energy per particle of the quasilattice being ~ 2%. 
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In the present paper we will consider the above optimal hard-sphere decoration of the 
Penrose tiling. Our treatment is based on the generalized effective liquid approximation!^ 
which has been previously applied to perfect hard spheres and hard disks crystals yielding 
very accurate results as compared to the simulation data. The quasilattice sums are calcu- 
lated in the 3D real space using a method which substantially improves convergence errors 
in comparison to previous 6D reciprocal lattice treatments. 

The paper is organized as follows. In section [TT] we briefly review the generation of the 
3D Penrose tiling which allows us to introduce the hard-sphere decoration. Section pi 
summarizes the generalized effective liquid approximation for the determination of the free 
energy of the quasicrystal. Our results are presented in section |TV] together with a discussion 
concerning the evaluation of quasilattice sums, while in the final section |V] we gather our 
conclusions. 



II. THE ICOSAHEDRAL QUASILATTICE 

In what follows, we consider the icosahedral quasilattice obtained by projecting a subset 
(to be specified below) of the 6D simple cubic lattice onto a 3-D hyperplaneH X n . The 
orientation of X n relative to the lattice is determined by requiring that the projected vectors 
{e-jj} of the 6D basis vectors {e J } (j = 1,2, ...6) coincide with the six vertex axes of the 
icosahedron, i.e.: 



p jk e k 



(2.1) 



where a summation over repeated labels is understood. The matrix representation of the 
projection operator P n is given by: 



/ 1 1 1 1 1 \ 

1 V5 1-1-1 1 
1 1 V5 1-1-1 
1-1 1 V5 1-1 
1-1-1 1 y/E 1 
VI 1-1-1 1 Vf> J 

and an elementary calculation leads to |e]| 2 = 1/2 (j = 1,2, ...6) and 

cos(e},ejj) = V5/5 (j = 2, 6) 



(2.2) 



cos(e,| J , e 



2+j a 2+k\ 



y/E/5 j-k = ±1, ±4 

-V5/5 j — k = ±2, ±3 



showing that the vectors {e-jj} may be identified with the six vertex directions of the icosa- 
hedron. 

We also consider the 3D hyperplane perpendicular to Xy, X x , obtained upon projection 
of the 6D basis vectors {e 7 } by the complementary projector P ± : 



e j _ pjkgk. pjk _ fijk _ pjk 



(2.3) 
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where 5 3 is the Kronecker delta. It can be readily shown that the projected vectors el 
may also be identified with the six vertex directions of the icosahedron, but permuted with 
respect to the projected vectors e^, i.e. ei • = — ejj • ef? (j ^ k). 

Both projections are dense in the 3D space but a quasilattice of finite density can be 
constructed by projecting onto X n only those points whose perpendicular space projection 
lies within a bounded region % known as the acceptance domaincl. To construct this bounded 
region we take the twenty distinct triplets {e^, ei, e^} of the projected vectors in X ± . Each 
triplet defines a rhombohedron of volume = |e^_ x ei • e*|. It is easily found that 
half of these rhombohedra are "large", i.e. = v^8sin(27r/5)/10 and half "small", i.e. 
Vijk — Vosin(47r/5)/10. The disjoint union of these twenty rhombohedra defines a closed 
convex region \ named the triacontahedron of volume v = \/8 [sin(27r/5) + sin(4-7r/5)]. The 
triacontahedron is therefore the projection onto X ± of the unit cell of the 6D cubic lattice. 

The selection of a subset of the 6D simple cubic lattice is accomplished by requiring that 
the orthogonal projections iri} of the lattice points {r 7 } of the 6D simple cubic lattice lie 
within the triacontahedronEa This construction yields a quasilattice in X^ which may also 
be regarded as a tiling of the 3D space by two kinds of rombohedrail. Using a one-parameter 
Gaussian approximation for the density peaks of the quasilattice, the one-particle density 
can be written as: 

PW=g) 3/2 E^(rl)e- Q(r - r " )2 (2.4) 

where the sum runs over the Bravais lattice vectors of the 6D simple cubic crystal, a is the 
inverse width of the Gaussians and the weight function iy(ri) is defined by: 

^( r± ^ I otherwise (2-5) 

It can be shown that the first three neighbour separations of the icosahedral quasilattice 
arei: r\ = (3 — 6y / 5/5)a 2 , r 2 = a and r| = (2 — 2y / 5/5)a 2 where a = y/2/2. Their 
average coordination numbers are 2/r 2 , 6 and 6, respectively, with r denoting the golden 
ratio r = (1 + y/E)/2. By locating a hard-sphere of diameter a — r\ at every vertex of the 
quasilattice, the packing fraction (the fraction of the total volume occupied by the spheres) 
is ~ 0.14, i.e. a packing fraction characteristic of a fluid phase. On the other hand, if we 
look for accomodating hard spheres of diameter a = a at every vertex, the short distance 
r\ does not allow this. But since the frequency of the nearest-neighbor separation is small, 
a better packing of hard spheres in the quasilattice can be obtained if one of the two vertex 
of each r!-bond is left vacant. These short bonds form closed rings of 10 links and chains 
of even or odd links. If we also require that two adjacent sites cannot be empty there are 
only two ways of placing the hard spheres on rings and chains. In the present investigation 
we have randomly located a hard sphere on a vertex of every ring and in an endpoint of 
each chain. This determines the accomodation of the remaining hard spheres, the effect of 
taking different initial localizations having a negligible effect on our results as the number of 
vertex of the quasilattice increases. It can be shown! that with this procedure, the volume 
occupied by the hard spheres increases leading to a hard-sphere packing fraction rj ~ 0.629, 
which is close to the random packing fraction (~ 0.64) and below the packing fraction of 
the crystal structures (~ 0.74 and ~ 0.68 for the fee and bec crystals, respectively). This is, 
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to our knowledge, the best icosahedral packing fraction of identical hard spheres reported 
so far. 



III. FREE ENERGY OF THE HARD-SPHERE ICOSAHEDRAL QUASILATTICE 

In recent years, the freezing of hard spheres into perfect crystals has been successfully de- 
scribed by several nonperturbative density functional theories. We here consider the stability 
of the hard-sphere icosahedral quasilattice described in |I| within one of such approaches, 
the generalized effective liquid approximation, which is now briefly summarized. 

The Helmholtz free energy F of a solid characterized by a one-particle density p(r) is a 
functional of p(r), denoted by F = F[p], which can be split as F[p] = F^[p] + F ex [p] where 

(3F id [p]= J rfrp(r)[ln{A 3 p(r)}-l] (3.1) 

is the ideal contribution with (5 = 1/ksT the inverse temperature and A the thermal de 
Broglie wavelength and 

(3F ex [p] = - J dvp{v) J dr'p(r') £ dX (1 - A) c(r, r'; [Xp]) (3.2) 



is the excess term. In ( |3.2|) c(r, r'; [Xp]) is the direct correlation function of the solid and 
X (0 < A < 1) is a parameter defining a linear path of integration in the space of density 
functions p\(r) = Xp(r) connecting a zero reference density to the one-particle density p(r) 
of the solid. The equilibrium solid density p(r), determined by functional differentiation, 
is the minimimum value of F[p] at constant average density. This variational calculation 
implies the direct correlation function of the solid which is the only unknown in (|3~T -j3~2) 
and hence some explicit approximations for i*ex[p] are required. 

Based on the similarity of the thermodynamic properties of the solid and fluid phases, 
the generalized effective liquid approximation first maps the excess free energy per particle 
of the solid onto that of some effective liquid, i.e.: 

i J drp(r) J dr'p(r') J dX (1 - A) c(r, r'; [Xp]) =p J dr J dX (1 - A) c(|r|; Xp) (3.3) 

where p is the density of the effective liquid, N — j drp(r) is the number of particles, and 



c(|r|; Xp) is the direct correlation function of the liquid. Equation ( |3.3|) is referred to as the 
thermodynamic mapping. 

In a second step, the generalized effective liquid approximation defines an structural 
mapping in which the direct correlation function of the solid is mapped onto that of a liquid. 
However, this mapping cannot be done directly because the direct correlation function of 
the liquid is translationally invariant while that of the solid is not. But taking into account 
that in ( |3.2j) the direct correlation function of the solid appears doubly weighted by the solid 
density , the difficulty is overcome by defining the structural mapping as: 

J drp(r) J dr'p(r')c(r,r';[p]) = J drp(v) J dr'p(r') c(\r - r'|; p[p]) (3.4) 
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With ( p.2H3.4| ) a self-consistent non-linear integral equation is obtained for the deter- 
mination of p[p] in terms of p(r) and the direct correlation function of the liquid. The 
complicated functional dependence p[p] can be simplified if p(r) is described in terms of a 
single order parameter a as in (|2.4|) in which case p becomes an ordinary function of a. 
The equilibrium solid density (i.e. a) is then determined by minimizing at constant average 
density the solid free energy with respect to the Gaussian width parameter a for a given 
crystal structure. 

As explained elsewhereS, the non-linear integral equation for the determination of p[p] 
can be further transformed into a system of two coupled nonlinear differential equations in 



z(X)-r/>(f)(X)) 



and z(X): 



z'(X) = 9(rj(\)) 



(3.5) 



(3.6) 



where 77(A) = np(X)a 3 /6 is the effective liquid packing fraction and a is the hard-sphere 
diameter. 

Using the one-parameter approximation fl2.4|) for the one-particle density of the quasi- 
lattice, $(r/(A)) is given by: 



1 N N oo 

^EEW^x) dRRc(R 
JV 1=1 j=l J0 



■,fj(X))S(R;a,ri 



(3.7) 



where m 



— r 3 ,, I and 



S(R; a, r. 





1/2 


a 







exp(-a(i? - r ij ) 2 /2) - exp(-a(R + r ij ) 2 /2) 



(3.8) 



In ( |3.5|) and (|3.6| ) the prime denotes the derivative with respect to the argument and 
ip(p)//3 is the excess free energy per particle of the fluid phase. For the latter we will use the 
Carnahan-Starling compressibility factor to obtain ip by thermodynamic integration of the 
equation of state while the Percus-Yevick equation is used for the structure of the fluid phase, 
i.e. the direct correlation function. Eqs. (|3.5|) and ( |3.6| ) have to be integrated numerically 
from A = to A = 1 with initial conditions rj(0) = z(0) = and the excess free energy per 
particle of the quasilattice ( |3.2| ) is finally determined as ^(^(l)). 



IV. RESULTS 



Before looking for a numerical solution of ( |3.5| - |3.6| ), we deal with a delicate point con- 
cerning the convergence of the quasilattice sums in ( |3.7|) . In order to emphasize it let us 
rewrite the rigth hand side of eq. (|3.7|) in the form: 



i N N 



jyEE^i) (4-1) 



=1 3=1 
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We first note that for a Bravais lattice ( f4.1|) reduces to, 

E zjMrj) (4-2) 
j 

where the sum runs over spherical shells of sites centered around the site at the origin, 
vj is the distance of shell j to the origin, and Zj is the number of sites at the jth-shell. 
In passing from (|4.1|) to (|4.2|) the translational symmetry of the Bravais lattice has been 
used. Moreover, since A(rj) decreases rapidly with distance, only a relatively small number 
of shells around the origin give a nonnegligible contribution to the sum leading to a rapid 
convergence of (|4.2|). 



But if the lattice does not have translational symmetry, the evaluation of ( 4.1]) becomes 
a delicate numerical problem if one looks for achieving a rapid convergence. For instance, if 
we consider a quasilattice of N sites the convergence of ( f4.1| ) (resulting from considering the 
N(N — 1)/2 different pairs i 7^ j) becomes so slow as N increases that it remains unreachable 
through usual computational efforts. This is also the case for a Bravais lattice as it may 



be tested by evaluating the sum ( |4.1| ) which, on the other hand, can be easily determined 
through 



As stated above, a possible way for dealing with the quasilattice sum! is to use the 6D 
reciprocal lattice. However, this procedure only provides a partial solution to the convergence 
of the sum since it is necessary to truncate the sum at some maximum value of the reciprocal 
lattice vector leading to an estimate error in the free energies of about 2%. 

We here propose an alternative solution to the convergence problem of ( |4.1| ) which pro- 



vides a substantial reduction of errors and computation time in the determination of the 
quasicrystal free energy. By starting with the 6D simple cubic lattice, we construct the 
quasilattice using the projection formalism described in [H[ Let M » 1 be the number of 
the lattice points generated in the quasilattice and draw an spherical surface containing al- 
most all the lattice points. Inside the sphere we construct a concentric sphere with iV >> 1 
(N < M) lattice points. 



By rewriting (|4.1|) as: 



1 N 
i=l 



5>( 



>:i 1 



(4.3) 



where we have separated the = contributions, we calculate the the i 7^ j-terms by first 
choosing a lattice point % (i = 1,2, ...N) inside the sphere and then summing consecutively 
over all neighbouring j lattice points until the relative error of the sum in brackets in ( |4.3| ) 
is less than a prefixed value. For lattice points % well inside the sphere, our procedure takes 
into account all the relevant "interactions" in ( f4.1| ). However, if the lattice point i lies near 
the surface of the sphere, our procedure overestimates the "interactions" in ( |4.1D because 
the sum in brackets contains points outside the sphere. We have found that these boundary 
errors can be reduced by increasing N (and therefore M) , in such a way that for M = 22000 
and iV = 10000, the estimate error of the free energy per particle of the quasilattice is about 
0.1%. Under such conditions, the computation time needed for evaluating the variational 
free energy (for each pair a-rj) is around 20' c.p.u. in a VAX 9000. 

In Fig. 1, the variational free energy per particle of the quasicrystal f3(f) = (3F/N — 
31n(A/a) + 1 is represented versus Y = (acr 2 /2) 1//2 for different packing fractions. We have 



7 



found minima for the quasicrystal free energy as a function of the Gaussian width parameter 
a for r] > 0.51, i.e. a stable or metastable quasicrystalline phase. Similarly to refB we have 
found that the quasicrystal turns out to be the more localized phase since the minima are 
always situated at greater F-values than the compact fee crystalline phase. Our results are 
gathered in Table | where the free energy and the Gaussian width parameter at the free 
energy mimimum of the quasicrystal are compared to those of the fee crystal. 

In Fig. 2 we represent the solid free energy per particle versus the packing fraction rj 
for the three solid phases (fee, bec and quasicrystal). We also include in the figure the fluid 
free energy obtained by thermodynamic integration of the Carnahan-Starling equation of 
state. It is seen that the quasicrystal is always metastable with respect to the remaining 
phases, the gaps of the free energy being somewhat less than those reported by M c Carley 
and Ashcroftfi 

V. CONCLUSIONS 

We have analyzed the stability of a hard-sphere quasicrystal obtained from a simple 
decoration of the 3D Penrose tiling which has been designed for optimizing the packing 
fraction. The stability has been analyzed using the generalized effective liquid approxima- 
tion. A simple method for evaluating the quasilattice sums in the 3D real space has been 
formulated. The method minimizes the boundary effects of finite quasilattices leading to a 
substantially better convergence than previous works. 

Our results show that the quasicrystal is metastable with respect to the crystalline and 
fluid phases. Such results agree with recent reported! calculations for a hard-sphere qua- 
sicrystal obtained from the modified weighted density approximation. Therefore, within 
these nonperturbative density functional theories entropy is insufficient to stabilize one- 
component quasicrystals. 

Since all known quasicrystals have complex metallic alloy phase structures it has been 
argued that for the stability of quasicrystals it is neccesary to have at least two class of atoms. 
The generalization of the one-component quasicrystal structure to an ordered two-component 
structure has been investigated by M c Carley and AshcroftEl who concluded that small changes 
in the diameter ratio of the two-component hard-sphere quasicrystal are not a stabilizing 
factor. It should be expected that energetic contributions resulting from considering more 
realistic interactions would propitiate stability. This possibility seems to be ruled out using 
the well-known perturbation schemes when applied to quasicrystals in view of the great free 
energy differences between the crystalline and quasicrystalline hard-sphere phases. Thus, 
the stability of quasicrystal structures within the modern density functional approaches is 
at present an open question. 
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FIGURES 



FIG. 1. Variational free energy per particle of the quasilattice (3(f) = (3F/N — 31n(A/<j) + 1 vs 
the Gaussian width parameter Y = (a<r 2 /2) 1 / 2 for packing fractions r/ = 0.51 — 0.60, by steps of 
0.01 (from bottom to top). 

FIG. 2. Free energy per particle (3(p = (3F/N — 31n(A/<r) + 1 vs the packing fraction r\ for the 
solid phases: fee (dotted line), bec (medium dashed line) and quasicrystal (dashed line), and for 
the fluid phase (continuous line). 
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TABLES 



TABLE I. Free energy f3(f) = (3F/N — 31n(A/cr) + 1 and Gaussian width parameter 
Y = (acr 2 /2) 1 / 2 at the free energy minimum for the quasicrystal (q) and the fee crystal (c) phases 
at different packing fractions 77. 
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P<t> q 


Y q 




Y c 


0.51 


6.15 


6.5 


5.28 


5.9 


0.52 


6.43 


7.0 


5.48 


6.5 


0.53 


6.73 


8.0 


5.68 


7.0 


0.54 


7.04 


9.0 


5.88 


7.6 


0.55 


7.37 


10.0 


6.09 


8.2 


0.56 


7.74 


11.0 


6.31 


8.9 


0.57 


8.14 


12.0 


6.53 


9.7 


0.58 


8.59 


13.5 


6.76 


10.6 


0.59 


9.12 


15.0 


7.01 


11.6 


0.60 


9.74 


16.5 


7.26 


12.6 
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